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ABSTRACT 


A  three-dimensional  (3D)  finite  difference  simulation  has  been  conducted 
of  surface  waves  propagating  across  a  laterally  heterogeneous  velocity  and 
density  model  of  the  Barents  Shelf.  The  3D  model  encompasses  a  region  of 
1620  by  810  km  to  a  depth  of  210  km.  The  calculations  were  designed  to  model 
0.02  to  0.05  Hz  surface  waves  from  an  isotropic  source  at  Novaya  Zemlya  in  a 
3D  heterogeneous  model. 

Extensive  focusing/defocusing  of  surface  waves  are  observed  for  the 
entire  bandwidth  from  0.02  to  0.05  Hz.  Focusing/defocusing  is  observed  to  be 
frequency  dependent  such  that  the  amplitude  patterns  of  0.02  and  0.05  Hz  are 
not  correlated.  A  hybrid  Fresnel-Kirchoff  technique  is  used  to  project  Rayleigh 
wave  amplitudes  to  teleseismic  distances.  About  25  percent  of  the  variance  of 
the  far-field  Rayleigh  wave  amplitude  scatter  from  Novaya  Zemlya  events  can 
be  explained  by  near-source  (within  800  km)  scattering. 

Multipathing  was  observed  to  occur  for  0.04  to  0.05  Hz  Rayleigh  waves, 
consistent  with  observations  of  Novaya  Zemlya  events  observed  at  NORSAR. 
Rayleigh-to-Love  conversion  is  observed  in  the  simulation  in  localized  regions 
along  structural  boundaries.  Such  scattered  Love  waves  should  be  observable 
at  narrow  azimuthal  ranges. 


1.  INTRODUCTION 


Surface  waves  from  Novaya  Zemlya  encounter  substantial  3D  lateral 
heterogeneity  as  they  propagate  in  the  first  several  hundred  kilometers.  The 
transition  from  thickened  continental  crust  to  shelf,  ocean,  and  shield  in  the 
vicinity  of  Novaya  Zemlya  is  particularly  complex.  The  complexity  of  the  region 
is  evident  in  the  LR  multipathing  observed  by  Levshin  and  Berteussen  (1979) 
for  Rayleigh  waves  arriving  at  NORSAR  from  Novaya  Zemlya  across  the  mixed 
Barents  Shelf  and  Fennoscandia  shield  path.  The  Barents  Shelf  to  the  west  of 
Novaya  Zemlya  has  deep  sedimentary  basins  that  are  responsible  for 
anomalously  slow  group  velocities  for  20  to  30  second  waves  as  observed  in 
northern  Fennoscandia  (Levshin  and  Berteussen,  1979;  and  Chan  and  Mitchell, 
1985).  Propagation  to  the  northwest  must  cross  oceanic  crust  with  thick 
sedimentary  cover  and  then  back  to  continental  crust  in  the  vicinity  of 
Spitsbergen  and  Franz  Josef  Land.  For  propagation  paths  to  Europe, 
Greenland,  and  North  America  the  surface  waves  must  propagate  obliquely 
across  various  portions  of  these  structures.  Zeng,  et  a!.,  (1989)  have 
demonstrated  that  the  focusing-defocusing  by  these  near-source  structures  and 
other  lateral  variations  in  the  Arctic  Ocean  produce  substantial  variations  in 
surface  wave  amplitude  across  North  America. 

In  previous  work  McLaughlin,  et  a!.,  (1991)  modeled  near  source 
scattering  of  surface  waves  by  lateral  variations  surrounding  the  Amchitka 
Island  Test  Site  in  the  Aleutian  islands.  It  was  demonstrated  that  both  lateral 
refraction  and  Rayleigh  to  Love  conversion  was  possible  in  the  presence  of 
shallow  heterogeneity  typical  of  a  subduction  zone.  In  this  paper  we  simulate 
the  propagation  of  surface  waves  through  structure  within  800  km  of  Novaya 
Zemlya  using  these  same  methods.  A  3D  velocity  and  density  model  is 
constructed  based  on  regional  geologic  models  from  the  literature.  Then  3D 
elasto-dynamic  finite-difference  calculations  are  used  to  simulate  the  complete 
seismic-wavefield  from  an  isotropic  source  in  the  3D  laterally  heterogeneous 
model.  The  simulations  are  then  analyzed  for  implications  on  the  propagation  of 
surface  waves  from  Novaya  Zemlya  across  the  Barents  shelf. 

We  find  that  laterally  heterogeneous  structure  in  the  Barents  Shelf  region 
west  of  Novaya  Zemlya  is  sufficient  to  introduce  far-field  focusing/defocusing, 
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multipathing,  and  Rayleigh  to  Love  conversion  in  the  0.02  to  0.05  Hz  bandwidth 
consistent  with  observations.  These  3D  finite  difference  simulations  provide 
insights  into  the  LR  scattering  mechanisms  and  provide  estimates  for  magnitude 
variations  that  can  be  introduced  by  strong  crustal  contrasts  in  the  vicinity  of  the 
source.  These  scattering  effects  are  significant  even  for  50  second  period 
Rayleigh  waves. 
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2.  OBSERVED  LR  PATH  EFFECTS 


LR  log-amplitudes  from  explosions  are  typically  scattered  with  standard 
deviations  in  the  range  of  0.2  to  0.3  log-units.  Even  when  systematic  source- 
receiver  path  corrections  are  taken  into  account,  momentimb,  and/or 
moment:log-yield  regressions  reveal  that  considerable  scatter  remains  in  the 
excitation  and  propagation  of  long-period  seismic  energy.  Both  the  uncertainty 
in  path  corrections  and  the  presence  of  tectonic  release  contribute  to  this 
scatter.  We  briefly  examine  data  from  the  northern  Novaya  Zemlya  test  site  to 
illustrate  the  nature  of  this  scatter.  First  we  examine  path  effects  alone  and  then 
the  evidence  for  tectonic  release. 

Figure  1  shows  86  station  corrections  with  estimated  standard  errors 
inferred  from  a  suite  of  19  northern  Novaya  Zemlya  explosions  by  Stevens  and 
McLaughlin  (1988).  The  explosions  are  assumed  to  be  isotropic  and  after 
spectral  Rayleigh  wave  amplitudes  are  corrected  for  path  dispersion  and 
attenuation  (see  Stevens,  1986),  station  corrections  are  derived  from  the  station 
log-moment  residuals  with  respect  to  the  network  average.  These  log-moment 
station  corrections  have  a  range  of  1.03  (factor  of  10.7)  and  a  standard 
deviation  (o)  of  0.218  (factor  of  1.65)  and  represent  76%  of  the  total  variance  in 
Rayleigh  wave  log-moment  scatter  seen  from  Novaya  Zemlya  in  the  0.02  to  0.05 
Hz  bandwidth. 

The  station  corrections  are  shown  plotted  in  Figure  1  as  a  function  of 
azimuth  from  0  to  360  degrees  as  well  as  an  expanded  plot  from  315  to  360 
degrees.  On  initial  inspection  the  station  corrections  appear  to  be  nearly 
random,  with  little  or  no  azimuthal  pattern.  However,  groups  of  stations  spaced 
within  a  few  degrees  of  each  other  in  azimuth  show  consistent  variation.  This 
deterministic  azimuthally  dependent  path  effect  is  clearly  apparent  for  stations 
in  the  northwest  quadrant  (315  to  360  degrees).  Within  the  estimated  error  bars 
the  station  corrections  exhibit  a  smooth  variation  in  amplitude  variation. 

To  quantify  this  correlation  between  adjacent  stations,  we  plot  a 
smoothed  variagram  and  the  inferred  correlation  function  in  Figure  2.  The 
variagram  is  a  plot  of  the  squared  differences  between  pairs  of  station 
corrections  plotted  against  the  separation  in  azimuth  between  the  stations.  The 
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Novaya  Zemlya  LR  Station  Corrections 


South 


Azimuth  (Deg) 


Figure  1 . 


Novaya  Zemlya  long-period  Rayleigh  (LR)  wave  station 
corrections  from  Stevens  and  McLaughlin  (1989).  Station 
corrections  (log-amplitude)  and  estimated  standard  errors  are 
plotted  at  a  function  of  station  azimuth  from  0  to  360  degrees 
(BOTTOM).  The  azimuth  range  from  315  to  360  degrees  is  shown 
at  a  larger  scale  (MIDDLE).  Amplitude  factors  are  plotted  in  a  polar 
format  at  the  azimuth  of  the  station  from  the  source  (TOP). 
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Figure  2.  The  square  of  the  each  3:ation  correction  pair  difference  is 
normalized  to  the  variance  and  plotted  versus  the  station  pair 
separation  distance  (BOTTOM).  The  inferred  correlation  function 
(TOP)  is  derived  from  the  smoothed  variagram  (BOTTOM). 


variation  between  stations  is  normalized  to  the  variance  of  the  distribution, 
Vjj  =  (Cj  -  Cj)2  /  a2,  where  q  is  the  i'th  station  correction.  If  the  individual  station 
corrections  are  random  and  independent  then  the  variation  should  average  to 
2.  If  the  station  corrections  are  correlated  then  the  average  variation  will  be  less 
than  2.  If  two  stations  are  close  together,  16]  -  0jl  <  10  degrees,  the  average 
variation  is  significantly  less  than  2.  suggesting  a  valid  correlation.  The  inferred 

correlation  function,  Cij(  50  y)  =  1  -  ^Vy,  inferred  from  the  smoothed  variagram 

is  shown  in  the  upper  panel  of  Figure  2.  For  pairs  of  stations  located  within  two 
degrees  in  azimuth  of  each  other,  the  station  corrections  are  30%  correlated. 
The  correlation  function  decays  for  increasing  distance  to  nearly  zero  for  station 
pairs  separated  in  azimuth  by  more  than  10  degrees  suggesting  that  they  are 
independent. 

We  therefore  conclude  that  30%  of  the  variance  (randomness)  in  the 
station  corrections  is  due  to  deterministic  path  effects  that  vary  in  azimuth  across 
the  network  with  characteristic  lengths  of  about  10  degrees.  Since  there  are 
many  gaps  in  network  coverage,  the  deterministic  pattern  is  aliased  for  most 
azimuths.  This  deterministic  pattern  represents  30%  of  the  total  observed 
variance  or  an  rms  variation  of  0.12  (factor  of  1.32)  in  log-moment  station 
corrections.  Because  deterministic  patterns  in  the  station  corrections  could  be 
due  to  non-isotropic  radiation  as  well  as  source-receiver  path  effects  we  need  to 
look  at  the  possibility  that  the  station  corrections  contain  tectonic  release 
information. 
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3.  TECTONIC  RELEASE 

Tectonic  release  at  the  Novaya  Zemlyatest  site  has  been  investigated  by 
Tucker,  et  al.  (1989,  1990),  Stevens  and  McLaughlin  (1988),  Burger,  et  al. 
(1986),  and  others.  Figure  3  shows  the  result  of  moment  tensor  inversion  of  the 
Rayleigh  wave  amplitudes  of  four  northern  Novaya  Zemlya  explosions.  Two  of 
these  explosions  have  a  very  small  amount  of  tectonic  release.  In  fact,  given  the 
amount  of  scatter  in  the  data,  the  tectonic  component  of  these  explosions  is  not 
significantly  different  from  zero.  Only  one  of  these  events  (1973255),  has  a 
significant  amount  of  tectonic  release.  These  moment  tensor  inversions  were 
performed  without  station  corrections,  and  as  is  evident  from  the  figures,  the 
amount  of  scatter  remaining  after  the  tectonic  contamination  is  removed  is 
substantial. 

Although  tectonic  release  is  present  in  some  events,  the  largest 
contributor  to  the  total  variance  in  log-moment  is  path  related.  From  the  analysis 
described  in  the  previous  section,  we  find  that  70  percent  of  the  station 
correction  variance  is  uncorrelated  for  station  separations  greater  than  10 
degrees,  and  station  corrections  account  for  76  percent  of  the  total  variance  of  a 
population  of  19  events. 
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f  =  0.12,  strike  =  6.882,  Mi  =  8.634E  >16  f  =  0.21,  strike  =  34.558,  Mi  =  4.817E  +  16 


Figure  3.  Inferred  radiation  patterns  from  four  of  the  best  recorded 
explosions  at  the  northern  Novaya  Zemlya  test  site.  No  station 
corrections  are  used  prior  to  inversion  for  tectonic  release.  No 
significant  non-isotropic  component  to  the  radiation  pattern  is 
statistically  significant  except  for  event  1973255,  when  station 
corrections  are  applied  prior  to  inversion  for  tectonic  release. 
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4.  THE  MODEL 

We  see  from  far-field  surface  wave  amplitude  data  that  there  is  evidence 
for  both  tectonic  release  and  source-receiver  path  effects.  In  order  to  examine 
propagation  effects  of  3D  structure  within  800  km  of  the  Novaya  Zemlya  test  site, 
we  have  constructed  a  3D  velocity  and  density  model  based  on  regional 
geologic  models  and  then  conducted  a  3D  finite  difference  simulation  of  surface 
wave  propagation. 

The  velocity  and  density  model  for  the  Novaya  Zemlya,  Barents  shelf, 
Kola  peninsula  region  was  based  on  data  from  Clarke  and  Rachlin  (1990) 
which  was  compiled  from  the  geologic  modeling  and  deep  seismic  sounding 
literature.  Figures  4b  and  4c  show  contours  of  crustal  thickness,  granitic  layer 
isopachs,  regions  of  missing  granitic  layer,  and  regions  of  reduced  sedimentary 
cover.  The  Barents  shelf  is  a  region  of  continental  crust  that  has  undergone 
considerable  extension  and  thinning.  The  basins  have  filled  in  places  with  as 
much  as  20  km  of  sediments.  The  Barents  Sea  is  bounded  to  the  east  by 
Novaya  Zemlya,  to  the  north  by  the  Svalbard  platform,  Spitsbergen,  and  Franz 
Josef  Land,  to  the  southwest  by  the  Kola  and  Kanin  Peninsulas,  and  to  the 
south  by  the  Pechora  Basin.  Novaya  Zemlya  is  considered  an  extension  of  the 
Ural  fold  belt  and  is  a  region  of  thick  crustal  thickness  (40-50  km).  The  Kola 
peninsula  is  an  old  shield  area  of  thick  (40  km)  crust  and  reduced  sedimentary 
cover.  Similarly,  the  Svalbard  platform  and  regions  around  Spitsbergen  and 
Franz  Josef  Land  to  the  north  are  continental  crust  (40  km).  The  Pechora  Basin 
between  the  Ural  fold  belt  and  the  Kanin  peninsula  has  a  reduced  crustal 
thickness  (<  35  km)  and  is  related  to  the  extensional  features  of  the  Barents 
shelf. 

From  this  general  geologic  model  a  3D  seismic  velocity  model  was 
derived  for  the  region.  Velocities  and  densities  for  the  sediments,  "granite  layer," 
"basalti  ayer,"  and  mantle  were  each  given  a  depth  dependence  consistent 
with  the  velocity  section  shown  in  Figure  4d  adapted  from  Clarke  and  Rachlin. 
The  3D  velocity  model  was  sampled  on  a  grid  at  5  km  intervals  in  depth  and  at 
50  km  intervals  in  north-south  and  east-west  directions.  For  the  purposes  of 
calculation,  the  model  was  discretized  as  6  by  6  by  6  km  cells.  Cross  sections  of 
the  velocity  model  are  shown  in  Figure  5.  The  results  of  this  model  generation 
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Figure  4a  Map  of  the  Barents  Sea  region.  The  3D  finite  difference  grid 
modeled  long-period  wave  propagation  in  a  rectangular  region 
1620  km  by  810  km  around  the  northern  Novaya  Zemlya  test  site. 
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Figure  4b. 


1 
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Map  of  crustal  thickness  (Moho  depth  contours  in  km)  of  the 
Barents  shelf  region  based  on  seismic  and  gravity  data.  Stippled 
areas  are  regions  of  missing  granitic  crust.  Hatched  areas  are 
regions  of  reduced  sedimentary  cover.  Adapted  from  Clarke  and 
Rachlin  (1990). 
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Figure  4c. 
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Map  of  granite-metamorphic  layer  (isopach  contours  in  km)  of  the 
Barents  shelf  region  based  on  seismic,  gravity,  and  magnetic  data. 
Hatched  areas  are  regions  of  missing  granite-metamorphic  layer 
inferred  from  magnetic  data.  Adapted  from  Clarke  and  Rachlin 
(1990). 
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Figure  4cj. 


Cross  section  across  the  Barents  Shelf  region  derived  from 
seismic  and  gravity  data.  Adapted  from  Clarke  and  Rachlin  (1990). 
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Shear  Velocity  (Azimuth  225) 


< -Southwest  (km) 


Shear  Velocity  (Azimuth  270) 
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Shear  Velocity  (Azimuth  315) 


<- Northwest  (km) 


Figure  5.  Shear  velocity  contours  on  three  sections,  northwest  (BOTTOM), 
west  (MIDDLE),  and  southwest  (TOP). 
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was  a  pastiche  of  18  different  crustal  and  upper  mantle  velocity  models 
arranged  on  a  grid  (see  Figure  6).  Phase  and  group  velocities  for  these  18 
different  crustal  models  are  shown  in  Figure  7  to  demonstrate  that  they  fall  into 
roughly  five  different  sets  of  dispersion  curves  corresponding  to  shield,  shield 
with  sediment  cover,  shallow  sedimentary  basin,  deep  sedimentary  basin,  and 
oceanic  models.  Figures  8a  and  8b  show  fundamental  Rayleigh  wave  group 
velocities  for  several  frequencies  contoured  on  the  surface  of  the  model  grid. 
The  low  velocity  basins  are  clear  features  at  all  frequencies  but  most  dominant 
at  the  higher  frequencies  (0.05  Hz). 
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Figure?.  Phase  velocities  (TOP)  and  group  velocities  (BOTTOM)  of  18 
individual  regions  of  the  3D  model  of  the  Novaya  Zemlya  -  Barents 
Shelf  area.  The  18  regions  cluster  into  5  distinct  groups. 
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Figure  8a.  Contours  of  Rayleigh  wave  group  velocity  for  the  3D  model 
illustrate  the  complexity  associated  with  the  Barents  Shelf 
sedimentary  basins. 


Group  Velocity  at  0.026  Hz  Group  Velocity  at  0,034  Hz 


Group  Velocity  at  0.042  Hz  Group  Velocity  at  0.050  Hz 


Shaded  gray-scale  contours  of  group  Rayleigh  wave  velocity  for 
the  3D  model.  Note  the  increase  in  complexity  and  contrast  in  the 
group  velocity  structure  for  the  higher  frequencies. 


5.  FINITE  DIFFERENCE  CALCULATIONS 


The  computational  grid  was  designed  to  model  20  to  50  second  surface 
waves  from  an  isotropic  explosion  in  a  laterally  heterogeneous  structure.  The  6 
km  cell  size  ensured  that  the  shortest  shear  wavelengths  of  interest  were  longer 
than  8  to  10  cells.  The  useful  model  as  indicated  in  Figure  4a  is  a  grid  of  270  by 
135  by  35  cells  (1620  km  north-south,  810  km  east-west,  and  210  km  vertical). 
Experience  gained  from  similar  calculations  (McLaughlin,  et  al.,  1991)  showed 
that  a  200  km  deep  grid  was  adequate  to  propagate  20  to  50  second  period 
surface  waves.  An  east-west  reflection  symmetry  was  assumed  for  the  velocity 
model  across  a  north-south  line  through  the  source  so  the  effective  useful 
model  was  1620  by  1620  by  210  km.  Outside  the  useful  model  to  the  north, 
south,  west,  and  at  depth  were  buffer  regions,  15  cells  thick.  The  cell  sizes  in 
these  buffer  regions  increase  at  a  compounding  rate  of  5%  outside  the  useful 
model  region  to  place  the  grid  boundaries  farther  from  the  useful  grid  than 
would  be  possible  with  a  grid  of  constant  cell  size.  Artificial  viscosity  is  used  to 
damp  reflections  from  the  regions  of  expanding  grid  and  reflections  from  the 
grid  boundaries  to  the  north,  south,  west,  and  at  depth.  Therefore,  the  total 
computational  grid  was  301  by  151  by  51  for  a  total  of  2.318  million  cells.  The 
elastodynamic  finite  difference  equations  were  integrated  with  a  time  step 
(cycle)  of  0.25  seconds  for  a  total  of  500  seconds  (total  of  4.6  billion  cell-cychs). 
The  fully  vectorized  computation  using  the  TRES  3D  finite  difference  code 
required  approximately  3  hours  of  CPU  time  on  a  256  million  word  CRAY  2. 

An  explosion  source  is  introduced  as  an  isotropic  moment  tensor  with  a 
time  function  proportional  to  te-®i  and  a  =  0.1  Hz.  This  is  a  minimum-phase 
source  time  function  with  a  corner  frequency  of  0.1  Hz  and  high-frequency  roll¬ 
off  proportional  to  o)-2.  In  this  manner  the  simulation  is  dominated  by  the  low- 
frequency  bandwidth  of  interest  and  snap-shots  (or  movies)  of  the  displacement 
field  (surface  waves)  can  be  visualized  without  post-processing.  Displacements 
were  saved  for  every  fourth  time  step  on  the  surface  of  the  grid  (z=0)  and  on 
vertical  planes  passing  through  the  source  in  the  north-south  (x=0)  and  east- 
west  (y=0)  directions.  The  analysis  of  these  displacement  fields  is  presented  in 
the  next  section. 
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6.  ANALYSIS 


The  displacement  seismograms  saved  on  the  surface  of  the  3D  grid  at  1 2 
km  intervals  serve  as  a  basis  for  our  analysis  of  the  surface  wave  propagation  in 
the  3D  grid.  Snapshots  of  the  displacement  field  on  the  surface  of  the  3D  finite 
difference  grid  are  shown  for  190,  290  and  350  seconds  in  Figures  9a-c.  The 
vertical  (Z)  and  radial  (R)  components  of  motion  show  development  of 
asymmetrical  Rayleigh  wave  propagation  while  the  transverse  (T)  component 
shows  developing  Love  wave  coda.  Rayleigh  to  Love  conversion  is  most 
prominent  along  the  edges  of  the  sedimentary  basins  where  the  Rayleigh 
waves  must  propagate  at  oblique  angles  to  the  structural  contrasts.  Focusing  by 
the  low-velocity  sedimentary  basins  is  clearly  evident  on  the  vertical  and  radial 
components  of  motion.  Furthermore,  the  dispersed  wavefronts  are  distorted  as 
they  refract  around  the  low  velocity  basins  and  are  delayed  across  the  basins. 
This  is  clearly  visible  in  the  wavefield  to  the  southwest  and  northwest  of  the 
source. 

In  order  to  more  fully  understand  these  wave  phenomena  we  have 
processed  the  displacement  seismograms  on  the  surface  of  the  grid  using 
standard  surface  wave  processing  algorithms.  Narrow  band  filters  are  applied 
at  selected  frequencies  and  the  peak  amplitude,  phase,  and  group  delay  of  the 
Rayleigh  wave  is  estimated.  In  the  following  analysis  we  present  the  results  for 
25,  30,  40,  and  50  second  vertical  component  Rayleigh  waves. 

Seismograms  along  sections  to  the  west  and  southwest  are  shown  in 
Figures  10a,b,c  and  11a,b,c.  The  vertical  and  radial  displacement  waveforms 
show  a  Pn-Pg  wavefield  followed  by  the  much  larger  dispersed  fundamental 
Rayleigh  wavetrain.  Vertical  and  radial  motion  are  fairly  coherent  with  each 
other  while  the  transverse  component  has  a  weak  Pn-Pg  coda  followed  by  an 
emerging  coda  after  the  Rayleigh  arrival.  The  westerly  section  has  lower 
transverse  levels  than  the  southwestern  section.  Coherent  propagating 
waveforms  can  be  seen  propagating  on  the  transverse  component  in  the 
southwesterly  direction.  Peak  amplitudes  of  the  transverse  component  waves 
propagating  to  the  southwest  exceed  1/2  of  the  radial  component  peak 
amplitudes. 
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nt  field  on  the  surface  of  the  3D  finite  difference  grid.  The 
I  show  the  development  of  the  asymmetrical  Rayleigh  wave 
developing  Love  wave  coda. 


Figure  9b.  Snapshots  (at  290  seconds)  of  the  displacement  field  on  the  surface  of  the  3D  finite  difference  grid.  The 
vertical  (Z)  and  radial  (R)  components  of  motion  show  the  development  of  the  asymmetrical  Rayleigh  wave 
propagation  while  the  transverse  (T)  shows  the  Rayleigh  to  Love  conversion  that  occurs  along  the  edges 
of  the  sedimentary  basins. 


Figure  9c.  Snapshots  (at  350  seconds)  of  the  displacement  field  on  the  surface  of  the  3D  finite  difference  grid.  The 
vertical  (Z)  and  radial  (R)  components  of  motion  show  focusing  by  the  low-velocity  sedimentary  basins. 
The  transverse  (T)  shows  Rayleigh  to  Love  conversion  along  the  edges  of  the  sedimentary  basins. 


Figure  10a.  Vertical  displacement  seismogram  section  to  the  west  (270  degrees).  Relative  peak  amplitude 
values  are  shown  to  the  right  of  each  trace. 
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Figure  10c.  Transverse  displacement  seismogram  section  to  the  west  (270  degrees).  Relative  peak  amplitude 
values  are  shown  to  the  right  of  each  trace. 
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Figure  11a.  Vertical  displacement  seismogram  section  to  the  southwest  (220  degrees).  Relative  peak  amplitude 
values  are  shown  to  the  right  of  each  trace. 
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Figure  11b.  Radial  displacement  seismogram  section  to  the  southwest  (220  degrees).  Relative  peak  amplitude 
values  are  shown  to  the  right  of  each  trace. 
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Figure  11c.  Transverse  displacement  seismogram  section  to  the  southwest  (220  degrees).  Relative  peak 
amplitude  values  are  shown  to  the  right  of  each  trace. 


In  order  to  get  a  better  picture  of  the  frequency  dependence  of 
amplitudes  in  the  grid,  narrow  band  filters  were  applied  to  the  vertical 
component  seismograms  and  the  peak  amplitudes  are  contoured  for  25,  30,  40, 
and  50  second  periods  in  Figures  12a,b,c,d.  The  amplitudes  have  been 
multiplied  by  Vr  in  order  to  make  a  partial  correction  for  geometrical  spreading. 
Note  that  asymmetries  are  present  at  all  periods.  The  amplitude  contours  of  the 
short  period  waves  (25  and  30  seconds)  outline  some  of  the  deep  sedimentary 
basins.  The  longer  period  amplitude  contours  clearly  show  how  the  amplitudes 
are  larger  to  the  west  than  to  the  north  or  south. 

The  same  processing  was  repeated  for  the  radial  and  transverse 
components  of  motion.  The  peak  amplitude  contours  at  25  and  40  second 
periods  are  shown  in  Figures  13a,b  and  Figures  14a,b.  The  transverse  motion 
is  minimum  to  the  west,  north,  and  south  with  several  ridges  of  high  amplitude 
extending  to  the  northwest  and  southwest.  These  high  amplitude  ridges  of 
radial  and  transverse  components  at  25  second  period  correlate  with  the  radial 
extension  of  the  edges  and  corners  of  sedimentary  basins.  The  zones  of 
Rayleigh-to-Love  conversion  at  25  and  40  second  periods  are  best  seen  in 
Figures  15a,b  where  the  T/R  ratio  is  contoured. 

In  addition  to  the  peak  amplitudes  at  selected  periods,  the  processing  of 
the  seismograms  can  yield  the  phase  and  group  delays  at  each  location  on  the 
surface  of  the  grid.  Group  delays  are  shown  for  25,  30,  40,  and  50  second 
periods  in  Figures  16a,b,c,d.  Wave  energy  travels  normal  to  group  delay 
contours,  and  this  analysis  serves  to  show  that  the  direction  of  propagation  is 
complicated  to  the  northwest  and  southwest.  Even  the  long  period  (50  second) 
group  delay  contours  show  the  effects  of  lateral  refraction.  At  the  shortest  period 
(25  seconds)  large  perturbations  in  the  wavefront  can  be  seen  along  the 
southern  boundary  of  the  Barents  Shelf.  The  slow  structures  to  the  west  retard 
the  waves  while  the  faster  structures  to  the  south-southwest  allow  a  faster  route 
for  the  wave  energy.  Multipathing  at  25  second  periods  becomes  evident  where 
contours  close  in  on  themselves  to  the  southwest  and  northwest. 

Similarly,  phase  delay  contours  can  be  constructed  across  the  grid 
surface  and  the  average  phase  velocity  from  source  to  receiver  can  be 
contoured.  Figure  17a,b  shows  these  contours  of  phase  delay  and  inferred 
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<-West  East-> 


Figure  12c.  Contours  of  40  second  period  peak  vertical  Rayleigh  wave 
amplitude  on  the  surface  of  the  3D  grid.  Note  the  asymmetries 
associated  with  the  sedimentary  basins. 


<-West  East-> 


Figure  13b.  Contours  of  40  second  period  peak  radial  amplitudes  with  1/^^ 
geometrical  spreading  correction  on  the  surface  of  the  3D  grid. 
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Figure  14a.  Contours  of  25  second  period  peak  transverse  amplitudes  with  1/Vr 
geometrical  spreading  correction  on  the  surface  of  the  3D  grid. 
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Figure  14b.  Contours  of  40  second  period  peak  transverse  amplitudes  with  1/Vr 
geometrical  spreading  correction  on  the  surface  of  the  3D  grid. 
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<-West  East-> 


Figure  15b.  Contours  of  the  ratio  of  Transverse/Radial  40  second  period 
motion  on  the  surface  of  the  3D  grid. 
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Figure  16a.  Group  delay  contours  for  25  second  period  vertical  Rayleigh 
waves  on  the  surface  of  the  3D  grid.  Contours  show  wave  energy 
refracting  around  the  Barents  Shelf  sedimentary  basins  to  the 
southwest  of  the  source.  Complexity  of  the  contours  at  distances 
greater  than  500  km  from  the  source  suggest  multipathing  at  25 
second  periods. 


Phase  Velocities,  NZ,  Vertical,  T=zC  sec 
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Figure  17b.  Inferred  phase  velocity  contours  for  25  second  period  vertical 
Rayleigh  waves  on  the  surface  of  the  3D  grid.  Phase  delay  has 
been  converted  to  phase  velocity.  The  southeastern  edge  of  the 
sedimentary  basin  to  the  southwest  of  tne  source  can  be  seen  as  a 
gradient  in  the  phase  velocities. 


r 


phase  velocity  for  25  second  period  vertical  component  Rayleigh  waves.  The 
phase  velocity  contours  are  lowest  in  the  basins  with  gradients  serving  to 
highlight  the  edges  of  the  low-to-high  velocity  structural  contacts. 

Broadband  (0.04  to  0.05  Hz)  frequency- wavenumber  (fk)  processing  for 
two  25  station  arrays  on  the  surface  of  the  grid  are  displayed  in  Figure  18a,  and 
b.  The  array  in  the  southwest  corner  of  the  grid  (Figure  18a)  shows  Rayleigh 
wave  energy  (vertical  and  radial  components)  arriving  on  azimuth  while  the 
transverse  component  motion  arrives  nearly  13  degrees  off  azimuth.  The 
second  array  (Figure  18b)  was  placed  to  the  south-southwest  of  the  source  and 
records  vertical  and  radial  motion  only  7  degrees  off  azimuth  while  the 
transverse  motion  arrives  12  degrees  off  azimuth.  Both  arrays  are  recording  the 
Rayleigh-to-Love  conversions  from  the  sedimentary  basins  southwest  of  the 
source. 

To  summarize,  we  see  focusing/defocusing  of  surface  waves  with 
periods  ranging  from  25  to  50  seconds.  Refraction  and  multipathing  is  greatest 
for  the  shortest  wave  lengths.  Rayleigh-to-Love  conversion  is  occurring  where 
Rayleigh  waves  travel  obliquely  along  boundaries  between  high  and  low 
velocity  regions.  Transverse  motion  for  40  and  25  second  periods  can  locally 
exceed  radial  motion  but  regions  of  Rayleigh-to-Love  conversion  more  typically 
exhibit  T/R  ratios  between  1/4  to  1/2. 
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Figure  18a.  Broadband  (0.04  -  0.05  Hz)  frequency-wavenumber  (f-k)  analysis 
of  radial  (upper  left),  vertical  (upper  right),  and  transverse  (lower 
left)  component  motion  from  a  25  station  array  in  the  southwest 
corner  of  the  grid.  The  source  is  at  an  azimuth  of  45  degrees  from 
the  array.  Radial  and  vertical  energy  is  arriving  from  an  azimuth  of 
42  degrees,  transverse  energy  is  arriving  from  an  azimuth  of  58 
degrees. 
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Figure  18b.  Broadband  (0.04  -  0.05  Hz)  frequency-wavenumber  (f-k)  analysis 
of  radial  (upper  left),  vertical  (upper  right),  and  transverse  (lower 
left)  component  motion  from  a  25  station  array.  The  source  is  at  an 
azimuth  of  27  degrees  from  the  array.  Radial  and  vertical  energy  is 
arriving  from  an  azimuth  of  34  degrees,  transverse  energy  is 
arriving  from  an  azimuth  of  39  degrees. 
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7.  TELESEISMIC  LR  AMPLITUDES 


In  order  to  propagate  the  Rayleigh  waves  out  of  the  finite  difference  grid 
to  teleseismic  distances  we  need  to  make  some  simplifying  assumptions  about 
the  propagation  outside  the  finite  difference  grid.  In  order  to  do  this,  we  have 
applied  the  Fresnel-Kirchoff  integral  method  as  described  in  McLaughlin,  et  al. 
(1991)  to  predict  the  effects  of  near-source  propagation  (within  800  km)  upon 
the  amplitudes  of  teleseismic  long-period  Rayleigh  waves.  The  vertical 
component  Rayleigh  waves  observed  on  the  surface  of  the  3D  finite  difference 
grid  are  propagated  to  teleseismic  distance  (5000  km)  taking  into  account  far- 
field  diffraction  effects.  The  far-field  seismogram  in  the  frequency  domain, 
S(f,~r'),  is  given  by 

S(f,r)  e''^l'’~^(hi+h2)ndl  (1) 


where  Z  is  a  closed  loop  on  the  surface  of  the  finite  difference  grid,  s(f,  r )  is  the 
Rayleigh  wavefield  at  location  ?  on  the  loop  Z,  7  is  the  receiver  location,  k  = 
co/c,  c  is  the  far-field  phase  velocity,  is  the  normal  vector  from  r  to  r ,  A2  is 
the  propagation  direction  of  the  Rayleigh  wavefield  at  location  ~r,  and  i1  is  the 
normal  of  Z  at  location  7. 

It  should  be  noted  that  for  any  teleseismic  path,  additional  heterogeneity 
will  further  focus/defocus  the  surface  waves.  The  hybrid  Fresnel-Kirchoff  finite 
difference  (FK-FD)  method  only  models  near-source  scattering  as  it  would  be 
seen  in  an  otherwise  homogeneous  structure  outside  the  finite  difference 
model. 

The  predicted  FK-FD  far-field  Rayleigh  wave  amplitudes  in  several 
bandwidths  are  shown  in  Figure  19a  and  b  assuming  a  phase  velocity  range  of 
3.5  to  3.7  km/s.  The  effect  of  selecting  the  phase  velocity  to  propagate  Rayleigh 
waves  away  from  the  finite  difference  grid  is  addressed  in  Figure  19c  and 
shows  that  a  significant  variation  in  assumed  phase  velocity  has  little  effect  on 
the  apparent  radiation  pattern. 
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Figure  19a.  Predicted  Fresnel-Kirchoff  amplitudes  for  teleseismic  Rayleigh 
waves  in  several  bandwidths.  Note  the  smaller  overall  variation  in 
amplitudes  at  lower  frequencies  (0.020-0.027  Hz)  versus  the 
higher  frequencies  (0.043-0.050  Hz).  The  average  amplitudes 
(TOP)  over  the  widest  bandwidth  (0.020-0.050  Hz)  have  the  least 
variation  as  a  function  of  azimuth  from  the  source. 
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Fresnel-Kirchoff  Log— Amplitudes 


Figure  19b  Perspective  mesh  of  the  predicted  Fresnel-Kirchoff  amplitudes  for 
teieseismic  Rayleigh  waves  as  a  function  azimuth  from  the  source 
(190-370  degrees)  and  frequency  (0.02-0.05  Hz).  Note  the  greater 
variation  in  the  apparent  radiation  pattern  at  high  frequencies. 
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Figure  19c.  Predicted  Fresnel-Kirchoff  amplitudes  for  0.02-0.05  Hz  teleseismic 
Rayleigh  waves  assuming  two  ranges  in  phase  velocity,  3. 5-3.7 
and  3.8-4.0  km/s.  The  apparent  radiation  pattern  is  insensitive  to 
the  assumed  far-field  phase  velocity  used  in  the  Fresnel-Kirchoff 
integral. 
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The  FK-FD  method  predicts  localized  anomalies  in  amplitude  of  as  much 
as  50%  over  azimuthal  ranges  of  several  degrees.  Averaged  over  all  azimuths 
the  rms  variation  ranges  from  0.05  to  0.1  rms  log-amplitude  depending  upon  the 
frequency  bandwidth  of  interest.  Note  the  smaller  overall  variation  in  amplitudes 
at  lower  frequencies  (0.020-0.027  Hz)  versus  the  higher  frequencies  (0.043- 
0.050  Hz).  Amplitudes  averaged  over  the  largest  bandwidth  (0.020  to  0.050  Hz) 
have  reduced  variation  compared  to  the  any  sub-bandwidth. 

The  high-  and  low-frequency  surface  waves  are  affected  differently  by  the 
laterally  varying  structure  and  are  focused/defocused  at  different  azimuths.  For 
source  size  estimation,  averaging  over  the  widest  possible  bandwidth  may  be 
desirable  to  reduce  near-source  focusing-defocusing  effects.  Furthermore, 
given  the  greater  variability  of  the  apparent  radiation  pattern  at  high  frequencies 
(Figure  19b)  it  may  be  desirable  to  apply  a  frequency  dependent  weighting  to 
Rayleigh  wave  amplitudes. 
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8.  CONCLUSIONS  AND  DISCUSSION 


A  3D  geo-tectonic  model  based  on  review  of  seismic,  gravity,  and 
magnetic  surveys  was  used  to  construct  a  3D  computational  velocity  and 
density  model  for  a  3D  finite  difference  simulation.  The  calculation  was 
designed  to  simulate  20  to  50  second  waves  from  an  isotropic  explosion  source 
in  a  3D  heterogeneous  structure  1620  by  1620  by  210  km.  Although  2D  finite 
difference  simulations  are  becoming  commonplace  in  seismological 
investigations,  3D  computations  are  just  beginning  to  be  useful.  Such 
calculations  can  be  considered  to  scale  with  size  and  bandwidth,  so  a  similar 
calculation  could  be  performed  in  the  0.4  to  1  Hz  bandwidth  for  a  structure  81  by 
81  by  21  km.  Care  must  be  exercised  in  choosing  the  target  bandwidth  and 
designing  the  grid.  However,  only  a  3D  finite  difference  calculation  is  capable  of 
simulating  the  complexity  of  wave  propagation  phenomena  that  was  exhibited 
in  this  calculation.  Refraction,  focusing-defocusing,  and  wave  conversion  were 
observed  to  result  as  waves  were  scattered  by  heterogeneities  with 
characteristic  sizes  comparable  to  the  wavelengths  of  interest. 

Focusing-defocusing  is  observed  from  0.02  to  0.05  Hz  and  shows 
considerable  dependence  upon  frequency.  While  the  0.05  Hz  surface  waves 
are  strongly  affected  by  the  low  velocity  basins,  the  0.02  Hz  waves  see 
compensating  effects  of  a  shallower  Moho  and  are  refracted  differently.  These 
compensating  effects  suggest  that  focusing-defocusing  should  not  be 
considered  a  broadband  effect  and  averaging  over  the  widest  possible 
bandwidth  is  desirable.  Furthermore,  since  focusing-defocusing  of  surface 
waves  is  most  variable  at  shorter  wavelengths,  a  frequency  dependent 
weighting  may  be  appropriate  for  source  size  estimation. 

Although  there  is  ample  opportunity  for  teleseismic  surface  waves  to  be 
further  refracted  outside  the  near-source  region  the  Fresnel-Kirchoff  far-field 
amplitudes  suggest  that  near-source  scattering  may  be  responsible  for  about 
0.1  rms  log-amplitude  variation  in  narrow  bandwidths  and  about  0.06  rms  log- 
amplitude  scatter  in  the  0.05  to  0.02  Hz  bandwidth.  This  compares  with  0.2  rms 
log-amplitude  scatter  observed  in  station  corrections  from  the  northern  Novaya 
Zemlya  test  site  or  about  25%  of  the  total  variance.  From  our  variagram  analysis 
of  the  station  corrections  we  concluded  that  about  30%  of  the  variance  in  the 
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station  effects  was  correlated  over  narrow  azimuthal  ranges.  This  is  consistent 
with  the  simulation  results  that  predict  rapid  variations  in  teleseismic  amplitude 
with  azimuth. 

It  should  also  be  noted  that  estimation  of  attenuation  by  spectral  slope 
can  be  disturbed  by  frequency  dependent  focusing-defocusing  of  teleseismic 
surface  waves.  For  example,  some  azimuths  would  appear  as  low  Q  where  the 
0.02  Hz  waves  are  focused  and  the  0.05  Hz  waves  are  defocused.  At  other 
azimuths  the  reverse  can  happen  with  the  result  that  the  path  appears  as  high 
Q.  Focusing-defocusing  effects  may  therefore  contribute  to  uncertainties  in 
attenuation  path  corrections  and  further  contribute  to  the  scatter  in  station 
corrections. 

Evidence  of  surface  wave  multipathing  was  evident  in  the  3D  simulation 
as  surface  waves  found  faster  paths  around  low  velocity  basins  in  the  model. 
This  supports  the  multipath  observations  by  Levshin  and  Berteussen  (1979)  of 
NORSAR  recordings  from  Novaya  Zemlya  explosions.  However,  the  simulations 
suggest  that  these  multipaths  are  refracted  arrivals  rather  than  reflections  from 
structural  features  as  suggested  by  Levshin  and  Berteussen.  No  surface  wave 
reflections  were  clearly  apparent  in  the  simulations  despite  sharp  structural 
discontinuities  in  the  velocity  model.  Refracted  and  converted  energy  clearly 
dominates  the  scattered  wavefield. 

Rayleigh  to  Love  wave  conversion  was  observed  in  the  simulation. 
Conversion  is  most  intense  in  regions  where  Rayleigh  waves  are  obliquely 
incident  upon  structural  boundaries.  These  regions  were  most  clearly  seen 
where  the  wavefront  curvature  was  greatest  across  a  fast-to-slow  structural 
boundary.  The  Love  waves  were  then  observed  to  propagate  independently  of 
the  Rayleigh  waves.  In  localized  regions  of  the  grid,  transverse  to  radial 
component  ratios  exceed  1/4  at  0.04  Hz  and  1/2  at  0.025  Hz.  Frequency 
wavenumber  processing  confirms  that  Love  wave  type  motion  is  propagating  in 
the  grid  away  from  the  source  region.  As  in  a  related  study  (McLaughlin,  et  al., 
1991)  of  propagation  in  a  3D  subduction  zone  structure,  the  simulations  predict 
that  teleseismic  Love  waves  should  be  observed  in  narrow  isolated  azimuthal 
ranges.  It  is  quite  possible  that  such  weak  Love  waves  observed  on  a  sparse 
network  could  be  interpreted  as  Love  waves  from  tectonic  release. 
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ATTN:  Dr  L.  Turnbull 
OSI/NED,  Room  5P0130 
Washington  DC  20505 

I  US  Arms  Control  and  Disarmament  Agency  1 

ATTN:  Dr  M.  Eimer 

Verification  and  Intelligence  Bureau,  Rm  5741 
Washington  DC  20451 

US  Arms  Control  and  Disarmament  Agency  1 

I  ATTN:  Mr  R.J.  Morrow 

320  21st  Street  NW 
Washington  DC  20451 

US  Geological  Survey  1 

ATTN;  Dr  T.  Hanks 

I  National  Earthquake  Research  Center 

345  Middlefield  Road 
Menlo  Park  CA  94025 
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OTHER  GOVERNMENT  AGENCIES  (continued) 

US  Geological  Survey  MS-913  1 

ATTN:  Dr  R.  Masse 

Global  Seismology  Branch 

Box  25046,  Stop  967 

Denver  Federal  Center 

Denver  CO  80225 


UNIVERSITIES 


Boston  College  1 

ATTN:  Dr  A.  Kafka 

Western  Observatory 

381  Concord  Road 

Weston  MA  02193 

California  Institute  of  Technology  1 

ATTN:  Dr  D.  Harkrider 

Division  of  Geological  and  Planetary  Sciences 
Seismological  Laboratory 
Pasadena  CA  91125 

Columbia  University  1 

ATTN:  Dr  L.  Sykes 

Lament -Doherty  Geological  Observatory 
Pal isades  NY  10964 

Cornell  University  1 

ATTN:  Dr  M.  Barazangi 

INSTOC 

3126  Snee  Hall 
Ithaca  NY  14853 

Harvard  University  I 

ATTN:  Dr  J.  Woodhouse 

Hoffman  Laboratory 

20  Oxford  Street 

Cambridge  MA  02138 

Massachusetts  Institute  of  Technology  3 

ATTN:  Dr  S.  Soloman,  Dr  N.  Toksoz,  and  Dr  T.  Jordan 
Earth  Resources  Laboratory 
42  Carleton  Street,  E34-440 
Cambridge  MA  02142 

Southern  Methodist  University  2 

ATTN:  Dr  E.  Herrin  and  Dr  B.  Stump 
Geophysical  Laboratory 
Dallas  TX  75275 
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UNIVERSITIES  (continued) 


State  University  of  New  York  at  Binghamton  1 

ATTN;  Dr  F.  Wu 

Department  of  Geological  Sciences 
Binghamton  NY  13902-6000 

St  Louis  University  2 

ATTN:  Dr  B.  Mitchell  and  Dr  R.  Herrmann 

Department  of  Earth  and  Atmospheric  Sciences 

3507  Laclede 

St.  Louis  MO  63156 

Pennsylvania  State  University,  The  1 

ATTN:  Dr  S.  Alexander 

Geosciences  Department 

537  Deike  Building 

University  Park  PA  16802 

University  of  Arizona  1 

ATTN:  Dr  T.  Wallace 

Department  of  Geosciences,  Bldg  77 

Tucson  AZ  85721 

University  of  California,  Berkeley  1 

ATTN:  Or  T.  McEvilly 

Department  of  Geology  and  Geophysics 

Berkeley  CA  94720 

University  of  California,  Los  Angeles  1 

ATTN:  Dr  L.  Knopoff 
405  Hilgard  Avenue 
Los  Angeles  CA  90024 

University  of  California,  San  Diego  1 

ATTN;  Dr  J.  Orcutt 

Scripps  Institute  of  Oceanography,  IGPP,  A-025 
La  Jolla  CA  92093 

University  of  California,  Santa  Cruz  1 

ATTN;  Dr  T.  Lay 

Earth  Sciences  Board,  Institute  of  Tectonics 
Santa  Cruz  CA  95064 

University  of  Colorado  1 

ATTN:  Dr  C.  Archambeau 
CIRES,  Campus  Box  449 
Boulder  CO  80309 

University  of  Ill inois  1 

ATTN:  Dr  S.  Grand 

Department  of  Geology 

1301  West  Green  Street 

Urbana  IL  61801 
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UNIVERSITIES  (continued) 


University  of  Southern  California  1 

ATTN;  Dr  K.  Aki 

Center  for  Earth  Sciences 

University  Park 

Los  Angeles  CA  90089-0741 

DEPARTMENT  OF  DEFENSE  CONTRACTORS 

Analytical  Sciences  Corporation,  The  1 

ATTN:  Document  Control  (Dr  Richard  Sailor) 

55  Walkers  Brook  Drive 
Reading  MA  01867 

Applied  Theory,  Inc.  1 

ATTN;  Dr  J,  Trul io 

930  South  La  Brea  Avenue,  Suite  2 

Los  Angeles  CA  90036 

Center  for  Seismic  Studies  2 

ATTN:  AFTAC  Representative 
1300  N.  17th  Street,  Suite  1450 
Arlington  VA  22209 

ENSCO,  Inc.  1 

ATTN:  Mr  John  R.  Stevenson 
PO  Box  1346 
Springfield  VA  22151 

ENSCO,  Inc.  1 

ATTN;  Dr  R.  Kemerait 
445  Pineda  Court 
Melbourne  FL  32940-7508 

Gould,  Inc.  1 

ATTN:  Mr  R.J.  Woodard 

Chesapeake  Instrument  Division 

6711  Baymeado  Drive 

Glen  Burnie  MD  21061 

Maxwell  Laboratories,  Inc.  1 

S-CUBED  Reston  Geophysics  Office 

Reston  International  Center 

ATTN:  Mr  J.  Murphy,  Suite  1212 

11800  Sunrise  Valley  Drive 

Reston  VA  22091 

Pacific  Sierra  Research  Corp.  1 

ATTN:  Mr  F.  Thomas 

12340  Santa  Monica  Boulevard 

Los  Angeles  CA  90025 


64 


DEPARTMENT  OF  DEFENSE  CONTRACTORS  (Continued) 


Rockwell  International  1 

ATTN:  B.  Tittmann 
1049  Camino  Dos  Rios 
Thousand  Oaks  CA  91360 

Rondout  Associates,  Inc.  1 

ATTN:  Dr  P.  Pomeroy 
PO  Box  224 

Stone  Ridge  NY  12484 

Science  Appl ications  International  Corporation  1 

ATTN:  Document  Control  (Dr  T.  Bache,  Jr.) 

10210  Campus  Point  Drive 
San  Diego  CA  92121 

Science  Horizons  2 

ATTN:  Dr  T.  Cherry  and  Dr  J.  Minster 
710  Encinitas  Blvd,  Suite  101 
Encinitas  CA  92024 

S-CUBED,  A  Division  of  Maxwell  1 

Laboratories,  Inc 
ATTN:  Dr  Keith  L.  McLaughlin 
PO  Box  1620 

La  Jolla  CA  92038-1620 

Sierra  Geophysics,  Inc.  2 

ATTN:  Dr  R.  Hart  and  Dr  G.  Mel  1  man 
11255  Kirkland  Way 
Kirkland  WA  98033 

SRI  International  1 

ATTN:  Dr  A.  Florence 
333  Ravensworth  Avenue 
Menlo  Park  CA  94025 

Teledyne  Industries,  Inc.  1 

Teledyne  Geotech  Alexandria  Laboratories 

ATTN:  Mr  W.  Rivers 

314  Montgomery  Street 

Alexandria  VA  22314-1581 

Woodward-Clyde  Consultants  1 

ATTN:  Dr  L.  Burdick 
PO  Box  93254 
Pasadena  CA  91109-3254 
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NON-US  RECIPIENTS 


Ministry  of  Defense  1 

ATTN:  Mr  Peter  Marshal  1 

Blacknest,  Brimpton 

Reading,  FG7-4RS 

United  Kingdom 

National  Defense  Research  Institute  I 

ATTN:  Dr  Ola  Dahl  man 
Stockholm  80 
Sweden 

NTNF/NORSAR  1 

ATTN:  Dr  Frode  Ringdal 

PO  Box  51 

N-2007  Kjeller 

Norway 

University  of  Cambridge,  Bullard  Laboratories  1 

ATTN;  Dr  Keith  Priestley 

Department  of  Earth  Sciences 

Madingley  Rise,  Madingley  Road 

Cambridge,  CB3  OEZ 

United  Kingdom 

OTHER  DISTRIBUTION 

To  b3  determined  by  the  project  office  9 


TOTAL  80 
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